The electron density: a fidelity witness for quantum computation

There is currently no combination of quantum hardware and algorithms that can provide an advantage over conventional calculations of molecules or materials. However, if or when such a point is reached, new strategies will be needed to verify predictions made using quantum devices. We propose that the electron density, obtained through experimental or computational means, can serve as a robust benchmark for validating the accuracy of quantum computation of chemistry. An initial exploration into topological features of electron densities, facilitated by quantum computation, is presented here as a proof of concept. Additionally, we examine the effects of constraining and symmetrizing measured one-particle reduced density matrices on noise-driven errors in the electron density distribution. We emphasize the potential benefits and future need for high-quality electron densities derived from diffraction experiments for validating classically intractable quantum computations of materials.


Introduction
In this study, we demonstrate the evaluation of electron densities of molecules using quantum computation.We also suggest that the electron density can be a potent validation tool of future quantum calculations, which may prove intractable to solve with conventional quantum chemistry.The study of electron densities is central to several elds of chemistry, physics, and materials science.The Hohenberg-Kohn theorem stipulates that the electron density uniquely denes ground state properties of electronic systems. 1 Through the Hellmann-Feynman theorem, 2 electron densities provide information on the forces acting within molecules. 3,4As one of the most information-rich observables in physical sciences, [5][6][7][8][9][10] the density lays the foundation for density functional theory (DFT), a formalism for predicting properties of many-electron systems. 11As experiments are the arbiter of truth, the buck oen stops with the electron density.
Electron densities can, importantly, be reconstructed from the renement of X-ray diffraction and scattering data, 9 using, e.g., multipolar models, [5][6][7][8]10 X-ray constrained wave functions, 12 or the maximum entropy method. 13One motivation for our work is that experimentally determined electron densities can prove useful for testing the accuracy of future quantum computations of materials; simulations of which may be unfeasible with conventional computers. Toay, it is oen preferable, cheaper, and faster to obtain information on electron distribution through conventional quantum mechanical calculations, e.g., by solving the Schrödinger equation at some level of approximation. 14To obtain highly accurate computational results (energies, densities, or other properties), ab initio quantum mechanical methods are the most reliable.Unfortunately, the electronic structure problem scales exponentially with system size.2 In 1980 and 1981, Benioff, 15 Manin 16 and Feynman 17 pointed out that quantum computers may offer a way forward, enabling larger and more reliable simulations of quantum systems.18 In recent years, several quantum devices have been developed and applied to compute energies of small molecules (see e.g., [19][20][21] ).No such calculation has, however, yet exceeded the accuracy or speed of a conventional quantum chemical calculation.So, just how far away are we from useful quantum computation of chemistry?To quantiably answer that question, it is necessary to perform benchmarking and validate the results of quantum computation in some manner.Thus far, such efforts have almost exclusively focused on a comparison against total energies obtained using state-of-the-art (and near exact) quantum chemistry.[22][23][24][25] However, whereas accurate energies for smaller molecules are available, such comparisons will not be possible if future quantum calculations of more complex systems are made possible.We emphasize at the onset that whether the crossing of such a technological threshold is made in the current noisy intermediate-scale quantum (NISQ) era or if it will require full fault tolerance quantum devices does not affect our main message or conclusions.Here, we suggest that the subtle variability of electron densities, accessible either computationally or experimentally, can act as a potent benchmark 26 for the quantum computation of materials.We draw the colloquial analogy to delity witnesses in the title, 27 by which we mean experimentally accessible observables whose values (here in terms of topological features) help to quantify the delity of a quantum calculation.
In what follows, we demonstrate calculations of electron densities (and their topological features) of molecular hydrogen (H 2 ) and lithium hydride (LiH) using quantum computers.These molecules are archetypical examples of fundamentally different chemical bonds (covalent and ionic).Simulation of quantum hardware is also employed to derive the electron density of larger molecules, the lithium dimer (Li 2 ) and hydrogen cyanide (HCN).The quantum volume 28 and noise level of the devices we use are insufficient to demonstrate any advantage over classical implementations (details of the hardware are provided in the Methods section and ESI †). 29However, they suffice for our goala rst proof-of-principle evaluation of electron densities using quantum devices.By comparing aspects of the electron density topology in these molecules, we showcase a different way to benchmark the quality of quantum hardware calculations of chemistry.[36][37]

Topological analysis of electron densities
Within the Born-Oppenheimer formalism, the average oneelectron density for a system of N electrons can be expressed as 38 rðrÞ ¼ N ð J * el ðr; r 2 ; .; r N ; RÞ J el ðr; r 2 ; .; r N ; RÞdr 2 ; .; dr N ; ( where r, R, and J el denote the electronic and nuclear coordinates and the electronic wave function, respectively.The electronic density is a real-valued scalar eld, lending itself to topological analysis and the extraction of topological features f c .One way to analyze such topology is to study critical points (CPs), i.e., locations in the density where the gradient vanishes, Vr = 0, or is undened. 39The critical points form a concise set of features that yield insight into the molecular structure.
To characterize and distinguish between critical points, we make use of the Hessian and its trace, the Laplacian (or curvature) of the density,  ; (2) For our analysis and discussion, it suffices to distinguish between critical points using so-called signatures, k, dened as the sum of the signs of the three eigenvalues of the Hessian. 40he critical points in molecules we will analyze are of two forms: k = −3 indicates local maxima, such as positions of nuclei, where all curvatures are negative.A non-nuclear attractor (NNA) is a rare example of a CP located at off-nuclei positions.We will discuss one NNA in Li 2 .In contrast, k = −1 corresponds to saddle points in the electron density, i.e., positions where the curvature in one direction is positive.The latter topological feature is commonly referred to as bond-critical points (BCPs) because they are oen (but not always) found between neighboring atoms that are chemically bound.A BCP is a point of lowest electron density along a path of highest electron density; the lowest point along a ridge connecting maxima.Because the density uniquely denes them, BCPs constitute suitable points of comparison between levels of theory and experiment.BCPs are furthermore useful for characterizing bonds in different ways. 41For example, the sign of the Laplacian (eqn (3)) at a BCP is an indication of local depletion (if positive) or concentration (if negative) of the electron density relative to its surroundings.A negative sign of V 2 r indicates a covalent bond, while a positive sign hints at an ionic (or closed-shell) type of interaction.
Topological analysis of electron densities (viz.QTAIM) also offers a way of dening atoms within molecules. 40Within QTAIM, atoms are identied with basins (Fig. 1b and c), nonoverlapping regions of space within which all gradient trajectories of the electron density terminate at the same local maximum (i.e., a k = −3 critical point).Note, therefore, that basins need not strictly be centered around nuclei but can be associated with NNAs.Neighboring basins are separated by zero-ux surfaces where Vr(r)$n = 0, and where n is the normal to the surface at r.We will use such basins to ascribe partial charge to atoms inside molecules, and we suggest such measures, along with r and V 2 r at BCPs, as examples of potent density-based quantum computational delity witnesses.

Quantum computation of electron densities
For quantum computation, the electronic structure problem is conveniently expressed within second quantization. 42Within this formalism, the electron density can be dened as where f * ps ðrÞ and f qs (r) correspond to the n spin orbitals and p and s denote spatial and spin indices, respectively.D pq is an entry in the one-particle reduced density matrix (1-RDM). 43In the following, we use p and q as spatial variables and s and s to index spin, in turn, denoted as a or b (i.e., s,s ˛{a,b}).The 1-RDM can be expressed as, where, ha † ps a qs i = hJja † ps a qs jJi.
In eqn ( 5), a † ps and a qs are the fermionic ladder (creation and annihilation) operators of electrons in spin orbital f qs (r) while jJi represents the wave function, dened as a linear combination of Slater determinants, In this work, 1-RDMs are constructed following measurements of parametrized quantum circuits (ansätze), representing the ground state of molecules of interest.Measurement of the 1-RDM scales as O(n 2 ), rendering a density-based delity witness approach computationally efficient.We will return to discuss the sensitivity of D pq to noise on quantum devices and how the outcome of calculations can be affected by enforcing physically motivated constraints on the 1-RDM.
Several algorithms 18,44,45 can be utilized to prepare the ground state solution, jJ GS i, of a molecule by encoding the fermionic chemistry problem 18,46,47 onto a quantum computer.We rely on the variational quantum eigensolver 48 (VQE) algorithm in this study, as it has been thoroughly used with the current generation of quantum devices.The VQE algorithm, outlined elsewhere, 21,22 leverages both quantum and classical computation to iteratively optimize a parameterized quantum circuit U(q).It minimizes the expectation value of the Hamiltonian, min q hJ 0 jU † ðqÞ ĤUðqÞjJ 0 i; where jJ 0 i is the initial reference state, usually the Hartree-Fock conguration.The ground state problem can be reduced to the electronic Hamiltonian, dened in second quantization as where h pq and V pqrs are the one-and two-electron integrals.We stress that our specic choices of algorithms and encoding procedures are not essential for the general case of calculating electron densities with a quantum computer.What is necessary is (a) identication of a state of interest (in our case, the ground state) by some quantum algorithm or quantum simulation, 49 and (b) reconstruction of the 1-RDM through measurements of ha † p a q i, following eqn (5).

Effects of noise
Noise is perhaps the single most dening characteristic of current quantum computing and the NISQ era. 29NISQ algorithms are oen hybrid in nature, dividing the computational load between both quantum and conventional hardware.In our case, the orbitals f ps (r) and f qs (r) of eqn (4) are precisely known functions represented on a conventional computer and assumed errorless.The quantum computer, in turn, stores the orbital occupations and phase.It is in the quantum computer where noise enters as uncertainty and measurement errors.As our interest lies with the 1-RDM, we will focus on the specic effects of noise on the elements D pq .By dening elementwise errors as 3 pq = D pq − D NF pq , where D NF pq is a noise-and error-free reference value, we can divide the effects of noise on the off-diagonal elements into two categories: diagonally symmetric (3 ij = 3 ji ) and asymmetric (3 ij s 3 ji ) errors.
To see how different kinds of noise may affect off-diagonal elements, we rst look at a Jordan-Wigner mapping of our fermionic creation and annihilation operators: We here use k to index our qubits, where each qubit k maps to a unique spin orbital, f ps .X and Y represent the corresponding Pauli gates, whereas Z k/ is the application of Pauli Z gates to all qubits k − 1, k − 2, ., 1.Note that eqn ( 9) and ( 10) are identical except for the sign (phase) of the Y-gate.In practice, this similarity means that the measurements of ha † ps a qs i and ha † qs a qs i will both perform the same measurements of Pauli strings k 5 Z k/ and Y k 5 Z k/ , with the phase introduced as a classical coefficient.Any measured diagonal asymmetries in the 1-RDM should, therefore, be solely due to insufficient sampling of the noisy state.Thus, provided independent measurements, and time-independent noise for a given set of circuit parameters, q, we can expect to sample the same noisy quantum state P(q) with every measurement.Here, we use independent measurements to mean single evaluations of the same quantum circuit that are not affected by previous evaluations.We phrase the expectation value in density matrix representation (not to be confused with a reduced density matrix), where P(q) is the parameterized density matrix.In this representation, the expectation value of an operator Ô is given as h Ôi = Tr[P Ô].In other words, provided that the above assumptions hold, any pair of real-valued off-diagonal elements must be equal, Tr[P(q)a † i a j ] = Tr[P(q)a † j a i ].While noise is an unavoidable part of contemporary quantum computing, symmetries and properties in the studied system can oen be used to gauge and combat errors.As such, it is of interest to study how the 1-RDM and derived properties thereof are affected by noise.And conversely, how enforcing known symmetries and properties of the 1-RDM mitigates the effect of noise.We have opted to study two important properties of the 1-RDM: its hermiticity, D pq ¼ D * qp ; and that its trace equals the number of electrons of the studied molecular system, 52 Tr(D) = n el . 50ecause the 1-RDM is Hermitian, we can, assuming realvalued entries, enforce the desired transpose (D pq = D qp ) symmetry on the measured 1-RDMs by averaging the corresponding off-diagonal elements as We also ensure particle conservation by normalizing the sum of all diagonal elements, Tr(D), to equal the total number of electrons, n el , Diagonal elements of the density matrix can be measured directly on quantum hardware, while off-diagonal elements require rotations to the measurement basis.While these basis rotations technically introduce more noise to the off-diagonal elements, these effects will be dominated by the noise incurred by preparing the (approximate) ground state wavefunction.This is because the circuit ansatz to prepare the ground state contains many more gates (including noisier 2qubit gates) compared to the measurement basis rotations.Thus, diagonal and off-diagonal elements are expected to be affected by similar degrees of noise and are subject to the same rescaling.We refer to the combined result of eqn ( 11) and ( 12) as a noise-mitigated 1-RDM.2][53] Thus, we emphasize that our aim is not to propose a novel mitigation strategy but rather to study its effect on a noisy 1-RDM and the resulting topology of the electron density.As our focus lies in evaluating topological qualities of the electron densities in the presence of noise and errors, we have opted for a conceptually straightforward mitigation strategy.As will be noted in the Results and discussion section, some noise effects on the measured number of electrons can already be avoided by using certain fermion-to-qubit encodings and qubit tapering.In particular, the use of parity encoding 46 allows us to preserve the parity of the number of particles of each spin species (a and b).For the case of lone a and b electrons, such as the singlet H 2 , the parity conservation coincides with particle conservation.Consequently, performing the parity transformation in such cases will effectively protect against any error in the particle or projected spin numbers.

Results and discussion
To evaluate the viability of the electron density's topological properties as a form of delity witness, we perform quantum chemistry calculations on both quantum and conventional computers.The quantum calculations rely on hardwareefficient ansätze 54 and minimal or small basis sets to expand the molecular wavefunction (see the Methods section).We use high-quality electron densities from conventional coupled cluster calculations at the CCSD/aug-cc-pVTZ level of theory as reference data.We note that reference densities can, in principle (up to a feasible limit), be obtained in many ways, including other costly ab initio methods, more affordable DFT functionals, diffraction data, or even wavefunction-tting experiments. 55entral to our work are analyses of the signicant contribution of noise on the quality of electron densities obtained with current quantum calculations.To that end, we not only compare with noiseless quantum simulations and near-exact conventional calculations, but we also implement and evaluate the efficacy of the above-described mitigation strategy for topological properties.To make our test set feasible on available quantum hardware, we rely on small basis sets (STO-3G for H 2 , LiH and Li 2 , and 6-31G for HCN).A frozen (1s) core approximation is used for atoms heavier than H, together with relatively small active spaces: H 2 (2,2), LiH (2,3), Li 2 (2,4), and HCN (4,4).The canonical molecular orbitals used to dene these active spaces are provided in the ESI (Fig. S5-S8 †).Two kinds of encodings are used to map the fermionic spin orbitals to qubits: parity encoding with two-qubit reduction for H 2 and LiH, and Jordan-Wigner encoding for Li 2 and HCN.While quantum results from H 2 and LiH are obtained from real quantum hardware, results for Li 2 and HCN are from simulations of quantum hardware that include a depolarizing noise model.Specics of the noise model and further computational details are available in the Methods section and the ESI.† Fig. 2 demonstrates that signicant errors can be present in both the diagonal and off-diagonal elements of the 1-RDMs.Similar conclusions have been drawn by Arute et al. 52 and Smart et al. 56 Note that an error along the diagonal is especially detrimental, as thisunless there is fortuitous error cancellationcan result in an incorrect (and unphysical) total number of electrons (Table 1).Even when the number of electrons is conserved by means of the encoding method, as in the case of H 2 , errors are still apparent along the diagonal of the 1-RDM.Table 1 summarizes how over-and underestimation of the true number of electrons is (a) present for larger molecules and (b) can be effectively corrected by rescaling the 1-RDM.
We want to reiterate that errors in the off-diagonal elements come in two avors: symmetric and asymmetric with respect to the diagonal.Note that the latter type of error only arises in the data from real quantum devices.The found asymmetry agrees well with the reduced number of measurement samples used in these calculations relative to the 10 6 samples applied in our quantum simulations.Physical hardware can also experience noise levels uctuating over time, breaking our previous assumption of time-invariant noise.
Moving to topological analysis of our electron densities, we look rst at atomic partial charges (Table 2).The partial charges attributed to each topological atom are evaluated by partitioning space into atomic QTAIM basins (Fig. 1b and c).Table 2 shows atomic partial charges derived from such topological atoms and how noise can affect the quantication of this important chemical concept.Also shown in Table 2 are the results following the application of our error mitigation strategy of rescaling and symmetrizing the 1-RDM.Whereas the latter approach generally improves results, its occasional failures can impart valuable lessons on how to improve chemically informed error mitigation.
One case of failure of our mitigation strategy lies with the NNA of Li 2 , a rare (but well-known) feature in the electron density where a local maximum is present between the two nuclei (Fig. 3).One explanation for the reduced accuracy of the NNA's partial charge is a substantial over-occupation due to noise of the valence p orbitals and the 2s* orbital (Fig. S4 †), all of which have but a small overlap with the NNA basin.In contrast, the valence 2s orbital exhibits a signicant reduction in population due to noise and overlaps prominently with the NNA basin.In other words, it is the incorrect relative lling due to noise of the 1p and 2s* orbitals over the 2s that effectively removes electrons from the inter-nuclear NNA basin in favor of the nuclear basins.Because the rescaling aspect of our error mitigation strategy (viz.eqn (12)) acts on all orbitals proportionally to their occupation, the 2s is scaled down to a larger extent, further exaggerating the error in this case.One example where our approach to error mitigation does little to correct topological features is HCN.The main reason for this is that the used active space is small and only includes the two bonding p and two anti-bonding p* orbitals.In contrast to our example of Li 2 , the p orbitals in HCN are polarized due to the difference in electronegativity of N and C. The p* orbitals are more pronounced on C and, consequently, have a larger overlap with the identied C basin.In contrast, the p orbitals have a larger N contribution and overlap slightly more with the N basin.Because the net effect of noise is an overpopulation of the antibonding p* orbitals (Fig. S4 †), the charge distribution is slightly skewed in favor of C (Table 2).A second reason for the negligible effect of our mitigation scheme in the case of HCN is that our calculations are based on sufficiently many, 10 6 , simulated measurements, which effectively removes asymmetric errors in the 1-RDM.The main advantage of our error mitigation approach in this example is ensuring the correct number of electrons (Table 1).Fig. 3 illustrates the effect of noise as electron density difference maps in planes across the different molecules, while Table 3 quanties noise (and mitigation) at selected critical points in the density.Because our calculations rely on a hybridquantum-classical algorithm, where orbitals are handled on a classical computer, electron densities are always constrained by the symmetries inherent to those orbitals.Such symmetries can also protect against noise along certain mirror planes.For example, Table 3 shows how the electron density at BCPs reported for HCN are all invariant to noise.This robustness is not general but a consequence of the small (4e, 4o) active space used.Because the p and p*orbitals have nodal planes along the bond axis, it is only in the Laplacian of r(r) that one can distinguish between noisy, mitigated, and noise-free results for these points within HCN.For Li 2 , the situation is reversed, and the actual worsening of results following adaptation of our error mitigation is attributed to the rescaling of the 2s occupation, which substantially affects CPs on the boundary of the NNA basin.
Our examples serve to illustrate the need to account for noise imbalance in future NISQ-computation of chemistry that targets high delity and accuracy.What we mean by noise imbalance is that some observables are more protected against errors than others.For example, in the context of electron densities, symmetries of orbitals in the chosen active space, the encoding, the ansatz, etc., can all affect relative noise when comparing

Methods
Calculations on H 2 were performed with a hardware-efficient circuit on the Särimner device of Chalmers, employing 50 000 samples.Calculations for LiH were performed on ibmq_quito utilizing a hardware-efficient two-local ansatz using 8192 samples. 61To reduce the effects of noise, readout mitigation was applied to all hardware calculations.Calculations of H 2 and LiH made use of parity encoding to reduce the required number of qubits and circuit depths.Calculations of Li 2 and HCN were performed on simulated devices using two consecutive layers of the ExcitationPreserving hardware-efficient ansatz available in Qiskit version 0.42.1. 61imulations utilized a depolarizing noise model and were performed over 10 6 samples.All calculations were additionally simulated without noise for comparison.
Our near-exact witness data was obtained with conventional quantum chemistry methods, using PySCF 62 at the CCSD/augcc-pVTZ level of theory.Topological analyses of electron densities were performed with Critic2. 63Additional computational and hardware details are provided in the ESI.†

Conclusions
Whereas tests of quantum computational accuracy are straightforward for small computational problems, for which conventional computation can be referenced, this will not always be so.Due to the rapid advances in quantum computational hardware, we may eventually face situations when it is not easy to validate whether a given problem has been solved to our satisfaction.The motivation behind this work is a proof-ofconcept for utilizing electron densities as delity witnesses, future-proof benchmarks, for the quality of the quantum computation of materials and molecules.To that end, we have demonstrated the rst topological analyses of electron densities inside molecules carried out with the assistance of quantum computers.Our work focuses on measuring one-particle reduced density matrices on real and simulated quantum hardware, how error mitigation can be applied to these entities, and how the resulting electron densities differ from near-exact reference calculations.The molecular systems studied herein are small enough to be treatable by current quantum hardware and small enough for comparisons with high-quality electron densities derived from conventional computations.Our examples are chosen to demonstrate the sometimes-detrimental effects of noise being imbalanced both spatially and with respect to 1-RDM entries.Therefore, we suggest that noise imbalance should be considered when designing or selecting error mitigation techniques, ansätze, encoding, and active spaces.Our work emphasizes the growing potential and need for high-quality (quantum) crystallography experiments, 64,65 i.e., experimental determination of electron densities near the diffraction limit.

Fig. 1
Fig.1(a) Five-step procedure for electron density-based witness estimation.A wavefunction is initially converged using a quantum algorithm.The corresponding one-particle reduced density matrix (1-RDM, D) is subsequently measured from the converged wavefunction.The electron density, r, is constructed from the measured 1-RDM.Topological analysis of the constructed electron density generates topological features f c , such as critical points.Finally, the quality of the quantum calculation can be evaluated by comparing topological features to a known (possibly experimentally determined) witness.(b) Gradient field lines of the electron density and atomic basins in LiH derived from quantum computation (QC).(c) Gradient field lines and atomic basins for the near-exact witness, derived from conventional quantum chemistry calculations at the CCSD/aug-cc-pVTZ level of theory.Green and blue circles indicate BCPs in data from the quantum calculation and witness, respectively.

Fig. 2
Fig. 2 Errors in 1-RDMs, comparing quantum devices (hardware and simulated) to exact (statevector) calculations.Color and intensity indicate the sign and magnitude of the errors, respectively.Hardware data (for H 2 and LiH) demonstrate higher asymmetry than simulations (of Li 2 , HCN) with depolarizing noise models.Only 1-RDM elements that correspond to the active space are shown.All bonds are aligned with the x-axis.

Fig. 3
Fig.3Electron density difference between noisy and noise-free (top), mitigated and noise-free (middle), and between mitigated and witness (bottom) results.The increased occupation of anti-bonding orbitals due to noise can be seen in all studied systems, resulting in a decreased electron density between nuclei.

Table 1
Number of electrons obtained from the trace of the 1-RDMs.For H 2 , the noisy calculation retains all electrons due to parity encoding and qubit tapering.Our calculations of larger molecules deviate noticeably from the correct number of electrons due to noise, effects that can be corrected by error mitigation.Estimated error bounds are given as a standard deviation based on the number of samples used in the measurements.The electron number for H 2 is not affected by the number of samples Chalmers Särimner device.b ibmq_quito device.c Simulation using a depolarizing noise model. a

Table 2
Atomic partial charges derived from a topological QTAIM analysis of the electron density.Rescaling and symmetrization of 1-RDM results generally improve computed partial charges, with the clear exceptions of the NNA of Li 2 a Chalmers Särimner device.b ibmq_quito device.c Simulation using a depolarizing noise model.
local properties at different points in space.We note that multiple approaches exist that could potentially mitigate noise imbalance (see, for instance, ref.57-60 and references therein).

Table 3
Values of the electron density r(r) and its Laplacian V 2 r(r) at critical points in the density of a selection of molecules, comparing noisefree, noisy, and error mitigated data against a near-exact witness a CCSD/aug-cc-pVTZ level of theory.